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Vertical  structure  of  bottom  Ekman  tidal  flows:  Observations, 
theory,  and  modeling  from  the  northern  Adriatic 
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[i]  From  September  2002  to  May  2003,  fifteen  bottom-mounted,  acoustic  Doppler 
current  profilers  measured  currents  of  the  northern  Adriatic  basin.  Tidal  fluctuations  at  all 
seven  of  the  major  Adriatic  frequencies  were  synthesized  from  a  response  tidal  analysis  of 
these  measurements.  Most  observed  tidal  current  ellipses  were  nearly  reversing,  but 
near  the  bottom,  tidal  current  ellipses  all  shortened  and  broadened,  semidiurnal  currents 
led  upper  water  column  currents,  and  diurnal  tidal  current  ellipse  orientations  rotated 
counterclockwise  toward  the  bottom.  Theoretical  solutions  for  a  tidally  forced,  bottom 
Ekman  layer  with  vertical  eddy  viscosity  of  the  form  Az  =  (3z  +  k  were  least  squares  fit  to 
the  observations.  Average  values  were  (3  =  3  •  10  4  m/s  and  k=  5  •  10  4  m2/s.  The  value 
of  k  was  important  in  matching  tidal  orientation  and  phase  changes,  and  a  nonzero  (3 
was  important  in  matching  tidal  amplitude  changes.  The  Navy  Coastal  Ocean  Model 
(NCOM)  and  the  Quoddy  model  were  also  compared  to  the  observations.  The  average 
RMS  errors  for  the  bottom  Ekman  layer  were  0.22  cm/s  for  the  best  fit  theory,  0.35  cm/s 
for  NCOM,  and  0.36  cm/s  for  Quoddy.  Az  structures  from  NCOM  and  Quoddy  show 
that  time  variation  in  Az  is  relatively  unimportant  for  Adriatic  tides.  The  bottom  shear 
stresses  from  theory  were  larger  in  magnitude  than  those  from  the  bottom  drag 
formulations  in  NCOM  and  Quoddy. 

Citation:  Book,  J.  W.,  P.  J.  Marlin,  1.  Janekovid,  M.  Kuzmin,  and  M.  Wimbush  (2009),  Vertical  structure  of  bottom  Ekman  tidal 
flows:  Observations,  theory,  and  modeling  from  the  northern  Adriatic,  J.  Geophys.  Res.,  114,  C01S06,  doi:10,1029/2008JC004736. 


1.  Introduction 

[2]  In  his  pioneering  paper  on  the  dynamics  of  tides  on 
the  north  Siberian  Shelf,  Sverdrup  [1927]  first  solved  the 
bottom  Ekman  theory  problem  with  tidal  forcing  and  found 
agreement  with  observed  tidal  structure  changes  near  the 
bottom.  Since  that  time,  the  topie  has  been  further  explored 
and  advanced  by  many  investigators.  A  nonexhaustive  list 
includes  the  mainly  theoretical  works  of  Fjeldstad  [1929], 
Prandle  [1982],  Soulsby  [1983],  and  Yasuda  [1987],  and  the 
works  of  Munk  et  ai  [  1 970],  Kundu  et  al.  [1981],  Maas  and 
van  Haren  [1987],  Lueck  and  Lu  [1997],  Oilman  and 
Wilson  [1998],  Tsimplis  [2000],  Werner  et  al  [2003a],  and 
Davies  et  al  [2004],  all  of  which  compared  theories  of 
vertical  current  variation  to  current  measurements. 

[3]  However,  despite  this  body  of  work,  quantitative 
evaluation  of  the  theoretical  vertical  changes  in  tidal  structure 
has  been  limited  by  laek  of  appropriate  measurements.  Both 
the  studies  of  Kundu  et  al  [1981]  and  Maas  and  van  Haren 
[1987]  verified  the  major  characteristics  of  tidal  Ekman 
theory,  but  certain  aspects  of  their  fits  were  less  than 

'Naval  Research  Laboratory,  Stcnnis  Space  Center,  Mississippi,  USA. 

2Ccntcr  for  Marine  and  Environmental  Research,  Ruder  BoSkovid 
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satisfactory  (e  g  ,  the  veering  of  ellipses  in  the  study  by 
Kundu  et  al.  [1981]  and  the  diurnal  constituent  comparisons 
in  the  study  by  Maas  and  van  Haren  [  1 987])  and  they  faced 
major  measurement  limitations  (i.e.,  vertical  current  mea¬ 
surements  of  limited  duration  from  anchored  ships  by  Kundu 
et  al.  [1981]  and  point  current  samples  at  limited  depths  by 
Maas  and  van  Haren  [  1 987]).  Both  Lueck  and  Lu  [  1 997]  and 
Ullman  and  Wilson  [1998]  found  good  agreement  between 
logarithmic  layer  theory  and  measured  velocity  profiles  over 
short  durations  but  conducted  individual  fits  in  time  with 
bottom  drag  coefficients  varying  over  tidal  cycles  at  sites 
where  Coriolis  accelerations  could  be  neglected  rather  than 
comparisons  to  tidal  Ekman  theory.  Recently,  acoustic 
Doppler  current  profiler  (ADCP)  measurements  have  been 
compared  to  tidal  Ekman  theory  by  Tsimplis  [2000]  and 
Davies  et  al  [2004],  but  neither  of  these  studies  had  mea¬ 
surements  within  10  m  of  the  seabed  and  therefore  missed  a 
region  of  strong  change  in  tidal  characteristics.  Of  the 
structure  that  was  observed,  Tsimplis  [2000]  explained  90% 
of  the  variance  of  along-strait  velocity  and  phase  using  a 
combination  of  frictional  theory  and  internal  mode  theory  to 
account  for  internal  tides.  The  3-D  model  with  a  quadratic 
friction  law  used  by  Davies  et  al.  [2004]  agreed  well  with 
available  measurements  for  semidiurnal  tides,  but  less  well 
for  diurnal  tides,  likely  because  of  inaccuracies  in  the 
measurements.  Werner  et  al.  [2003a]  compared  tidal  velocity 
observations  at  one  location  to  1-D  models  using  either  a 
combination  of  linear  and  constant  eddy  viscosities  or  using  a 
Mellor-Yamada  level  2.5  turbulence  closure  scheme  [Mellor 
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Figure  1 .  Bathymetry  of  the  north  Adriatic  and  tide  ellipses  from  vertically  averaged  currents.  M2  ellipses 
are  drawn  in  magenta,  K!  ellipses  are  drawn  in  red,  S2  ellipses  are  drawn  in  yellow,  Oj  ellipses  are  drawn  in 
blue,  Pi  ellipses  are  drawn  in  green,  N2  ellipses  are  drawn  in  cyan,  and  K2  ellipses  are  drawn  in  black.  The 
velocity  scale  is  given  in  the  bottom  left  comer.  Place  names  used  in  this  paper  are  labeled,  with  cities 
indicated  by  black  dots.  Moorings  are  located  at  the  centers  of  the  ellipses  and  named  according  to  the  label 
given  to  eaeh  ellipse  set. 


and  Yamada ,  1982].  They  found  good  agreement  between 
the  data  and  both  types  of  models  during  times  of  low 
stratification  for  M2  tides,  the  dominate  constituent  of  the 
region  and  the  only  one  they  utilized. 

[4]  Recent  mooring  measurements  from  the  northern 
Adriatic  present  a  new  opportunity  to  evaluate  tidal  variation 
in  the  bottom  boundary  layer.  Velocities  were  measured  by 
ADCPs  at  15  different  locations  spread  throughout  the  basin 
for  more  than  six  months  over  the  winter  [Book  et  al. ,  2007b] 
together  with  pressure  measurements  from  wave/tide  gauges 
(WTGs)  at  most  sites.  Depth  cell  sizes  were  1  m  or  smaller, 
measurements  covered  the  entire  water  column  with  the 
exception  of  the  bottom  blanking  and  surface  contamination 
zones,  and  measurements  at  all  sites  were  made  within  3  m  of 
the  bottom  or  closer.  In  addition,  results  from  a  3-D  finite- 
element  model  of  the  Adriatic  dedicated  to  tides  [Janekovic 
and  Kuzmic ,  2005]  and  results  from  a  3-D  finite-difTerenee 
model  of  the  Adriatic  with  tides  [Martin  et  al. ,  2006]  are 
available  for  comparison  with  these  measurements  and  with 
theory. 

[5]  The  Adriatic  Sea  is  an  arm  of  the  Mediterranean  Sea. 
It  may  be  represented  roughly  as  an  800-km-long,  150-km- 
wide  channel,  oriented  southeast-northwest,  open  at  the 
southeast  end  (Strait  of  Otranto),  with  the  bottom  sloping 
upward  toward  the  closed  northwest  end.  The  northern 


Adriatic  (defined  here  to  oeeupy  the  region  northwest  of 
Aneona  and  Zadar)  is  the  final  200  km  of  this  “channel,” 
where  depths  gradually  slope  from  70  m  in  the  southeast  to 
less  than  1 0  m  in  the  northwest  (Figure  1 ).  The  southwest 
(Italian)  side  of  the  sea  is  characterized  by  a  mild  bathymetry 
slope,  the  Po  River  Delta  and  associated  Po  River  plume, 
and  a  boundary  current,  the  Western  Adriatic  Current,  with 
a  typical  strength  of  10  em/s.  In  contrast,  the  northeast 
(Croatian)  side  of  the  sea  is  characterized  by  nearly  vertical 
dropoffs  at  the  eoast,  numerous  deep  bays  and  channels,  and 
diffiise  and  varying  currents.  Most  of  these  bays  and  channels 
on  the  northeast  side  are  nearly  isolated  from  the  main 
Adriatic,  but  Kvamer  Bay  opens  up  to  the  Adriatic  through 
a  30-km-wide  passage.  Figure  1  shows  some  of  the  main 
features  of  the  northern  Adriatic. 

[6]  There  has  been  considerable  theoretical  and  praetieal 
research  on  the  tides  of  the  Adriatic  (see  Cushman-Roisin  et 
al  [2001,  ehapter  7]  for  a  review  of  work  prior  to  2001), 
but,  as  in  many  coastal  areas,  direet  measurements  of  tidal 
currents  have  been  limited  by  technological  and  fishing 
pressure  restrictions.  Malacic  et  al.  [2000]  extended  earlier 
semidiurnal  applications  of  the  theory  of  Taylor  [1921]  by 
Hendershott  and  Speranza  [1971]  and  Mosetti  [1986] 
with  a  general  theory  of  gravity  and  topographic  waves  to 
explain  the  dynamics  of  both  the  semidiurnal  and  diurnal 
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Table  1.  Mooring  Positions  and  Depths 


Mooring 

Lalitudc 

Longitude 

Depth  (m) 

SS2 

43.835  l°N 

I3.3066°E 

25 

SS4 

43.8836°N 

1 3.3667°E 

46 

SS5 

43.9307°N 

13  4261 °E 

57 

SS6 

43.9956°N 

I3,5044°E 

66 

SS8 

44  2567°N 

I3,9053°E 

65 

SS9 

44  4102°N 

14  1748°E 

59 

SSI0 

44  481 2°N 

I4,2904°E 

51 

CP2 

44  46 1 0°N 

I2,855I°E 

42 

CP3 

44,5402°N 

I3.1245°E 

42 

KB1 

44,7507°N 

I4,02I3°E 

48 

VR1 

45.3I39°N 

1 2.5081  °E 

17 

VR2 

45.2789°N 

I2,6370°E 

25 

VR4 

45.1878°N 

1 3.0281  °E 

33 

VR5 

45. 1 249°N 

13  2837°E 

35 

VR6 

45.0581°N 

13  5360°F 

33 

tides  of  the  northern  Adriatic.  Recent  work  [Cushman- 
Roisin  and  Naimie ,  2002;  Janekovic  et  al. ,  2003;  Janekovic 
and  M.  Kuzmic,  2005;  Martin  et  al. ,  2006]  has  focused  on 
using  3-D,  high-resolution,  numerical  models  with  realistic 
topography  and  Mellor-Yamada  turbulence  closure  schemes 
to  simulate  the  Adriatic  tides,  and  then  to  validate  the 
components  of  these  simulations  with  available  observations. 
In  addition  to  their  focus  on  depth-averaged  tides,  Janekovic 
and  Kuzmic  [2005]  also  compared  simulated  tidal  vertical 
structure  with  ADCP  observations  at  one  station  and  found 
good  agreement.  Differences  were  attributed  to  imperfectly 
represented  vertical  mixing.  The  interaction  of  stratification 
and  tides  in  the  Adriatic,  with  implications  for  the  vertical 
structure  of  tidal  currents,  is  also  presently  a  topic  of 
investigation,  with  Chavanne  et  al  [2007]  citing  Po  River 
stratification  as  a  possible  explanation  for  model/data  dis¬ 
crepancies  near  the  Italian  coast  and  Mihanovic  et  al.  [2006] 
publishing  a  dedicated  study  of  internal  tides  in  the  Adriatic. 

[7]  This  paper  reexamines  time-dependent,  tidally  forced, 
bottom  Ekman  structure  using  the  data  from  15  ADCPs  and 
results  from  the  Navy  Coastal  Ocean  Model  (NCOM)  and 
the  Quoddy  ocean  model.  Section  2  describes  the  measure¬ 
ments,  section  3  describes  the  models,  section  4  describes 
the  bottom  Ekman  tidal  theory  used  in  this  paper,  and 
section  5  compares  vertical  tide  structures  from  observations, 
theory,  and  models.  Finally,  discussions  and  conclusions  are 
presented  in  sections  6  and  7. 

2.  Measurements 

[s]  From  September  2002  to  May  2003,  an  array  of  RD 
Instruments  (RD1)  Workhorse  Sentinel  broadband  ADCPs 
was  deployed  in  the  northern  Adriatic  as  part  of  a  Joint 
Research  Project  (JRP)  between  the  U.S.  Naval  Research 
Laboratory  (NRL)  and  the  NATO  Undersea  Research  Centre 
(NURC)  The  JRP  moorings  consisted  of  14  trawl-resistant 
bottom-mounted  ADCPs  [Perkms  et  al. ,  2000]  distributed 
along  portions  of  4  mooring  sections.  An  additional  upward 
looking  ADCP  was  mounted  near  the  base  of  a  meteorolog¬ 
ical  tower  as  described  by  Cavaleri  [2000].  These  mooring 
positions  are  shown  in  Figure  1  and  listed  in  Table  1  with 
their  mean  sea  level  depths.  The  full  mooring  sections  were 
populated  by  both  the  JRP  moorings  and  moorings  from 
several  international  partners  collaborating  on  the  study  of 
the  northern  Adriatic  [Lee  et  al. ,  2005].  In  addition  to  the 
ADCP  measurements  of  currents  throughout  the  water 


column,  bottom  pressure  (by  ADCP  or  wave/tidc  gauge) 
was  also  measured  at  each  site.  Book  et  al  [2007a]  provide 
further  details  of  the  mooring  instrumentation  and  Book  et 
al.  [2007b]  show  monthly  mean  and  storm-driven  currents 
observed  at  these  sites. 

2.1.  Tidal  Analysis 

[9]  The  JRP  ADCPs  were  set  to  measure  the  currents 
using  62-90  s  bursts  of  1  Hz  pings  every  1 5  min,  except  for 
the  first  half  of  the  VR1  deployment,  which  used  16  min 
bursts  of  0.5  Hz  pings  every  hour  Quality  control  steps  to 
exclude  bad  data  included  an  objectively  determined  velocity 
error  cutoff  (velocity  errors  estimated  from  independent 
measures  of  vertical  velocity),  exclusion  of  ensembles  with 
more  than  40%  (20%  for  surface  measurements)  of  the  data 
marked  bad  by  internal  RD1  checks,  and  additional  tests 
described  by  Book  et  al.  [2007a].  The  surface  echo  interfer¬ 
ence  zone  was  truncated  above  a  time-varying  level  deter¬ 
mined  from  a  time  series  of  sea  surface  height  constructed 
from  the  pressure  and  acoustic-backscatter-intcnsity  mea¬ 
surements.  Linear  compass  drifts  (less  than  4°)  in  some 
records  were  verified  to  be  false  trends  by  tidal  analysis 
and  corrected  by  small,  gradual  rotations  of  current  vectors 
Despite  the  lack  of  any  physical  evidence  of  instrument 
malfunction,  the  orientation  disagreement  between  the 
observed  tidal  ellipses  at  station  VR5  with  neighboring 
observed  tidal  ellipses  and  three  independent  modeling 
simulations  strongly  suggests  a  compass  error  is  present  in 
the  VR5  data.  Therefore,  the  currents  at  site  VR5  were  rotated 
28°  clockwise  to  align  with  modeled  strong-constraint, 
variational,  data  assimilation  predictions. 

[10]  Tidal  analyses  of  the  ADCP  data  were  done  individ¬ 
ually  for  all  depth  cells  using  the  Response  Method  of  Munk 
and  Cartwright  [1966]  on  the  15-min,  current,  ensemble 
time  series  (time  values  assigned  to  the  center  time  of  the 
measurement  bursts).  Gaps  were  introduced  into  the  ADCP 
records  by  the  quality  control  steps  summarized  in  the 
previous  paragraph,  so  a  method  was  used  for  2-D  interpo¬ 
lation  and  extrapolation  to  replace  missing  values  in  the 
ADCP  records.  Sensitivity  tests  show  that  the  ncar-surfacc 
tidal  solutions  (primarily  only  the  ADCP  depth  cell  nearest 
the  surface)  were  influenced  by  interpolation  and  extrapo¬ 
lation  method  choices,  but  solutions  at  all  other  depths  were 
very  insensitive  to  this  as  they  had  few  gaps. 

[11]  After  the  missing  values  were  filled  in,  the  Response 
Method  was  used  to  solve  for  estimates  of  pure  tidal  time 
series  and  for  approximations  to  the  harmonic  tidal  coef¬ 
ficients  for  the  Oj,  Pj,  Kj,  N2,  M2,  S2,  and  K2  constituents. 
Tides  were  not  calculated  from  surface  ADCP  depth  cells 
with  50%  or  more  data  marked  bad  by  the  quality  control 
procedures.  Slightly  different  procedures,  described  by  Book 
et  al.  [2007a],  were  used  to  calculate  tidal  currents  at  site  VR 1 
to  account  for  the  change  in  ADCP  settings  midway  through 
the  deployment. 

[12]  For  each  depth  level,  error  estimates  were  obtained 
by  creating  an  ensemble  of  40  normally  distributed,  random 
noise  time  scries,  each  having  variance  equal  to  the  ADCP 
measurement  error  variance.  Then,  the  velocity  time  series  of 
the  tide  and  of  the  tide  residual  were  obtained  by  response 
analysis  The  residual  was  2-h  low-pass  filtered  and  then 
added  to  the  tide  and  to  one  of  the  random  time  scries.  The 
procedure  was  repeated  for  each  random  time  scries  in  the 
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Figure  2.  Vertical  structure  of  M2  tidal  currents  for  all  sites  except  KB1  and  VR1  2.  Conventions  for 
tidal  current  parameters  are  those  of  Foreman  [1978]. 


ensemble.  95%  confidence  limits  for  the  tidal  ellipse  param¬ 
eters  were  assigned  to  1 .96  times  the  value  of  the  RMS  errors 
of  the  parameters  from  the  ensemble.  The  RMS  errors  were 
all  smaller  than  the  parametric  bootstrap  error  estimates  fora 
standard  harmonic  analysis  of  the  data  using  a  colored 
bivariate  noise  model  [Pawlowicz  et  al,  2002]. 

2.2.  Tidal  Observations 

[n]  f  igure  1  shows  the  tidal  ellipses  of  the  vertically 
averaged  currents.  The  eccentricity  of  the  ellipses  is  high 
(i.e.,  nearly  reversing  tidal  currents)  at  all  sites  except  for 
KB1  in  Kvamer  Bay  and  VR1  and  VR2  in  the  northwest 
comer.  Except  for  K.B1,  the  ellipse  major  axes  are  approxi¬ 
mately  aligned  with  the  Adriatic  axis.  There  is  an  increase  in 
tidal  currents  at  the  Istrian  coast;  site  VR6  had  the  strongest 
tidal  current  with  an  M2  semi  major  axis  of  10  cm/s.  M2  tidal 
currents  are  strongest  at  all  sites  (average  semi  major 
axis  7  cm/s),  but  S2  and  tidal  currents  also  play  prominent 
roles  (average  semi  major  axes  4  cm/s  and  3  cm/s,  respec¬ 
tively).  Oi,  Pi,  N2,  and  K2  are  all  much  weaker,  with  average 
semi  major  axes  of  1  cm/s. 

[14]  The  vertical  structures  of  the  tidal  currents  for  M2 
and  K,  are  shown  in  Figures  2  and  3,  respectively.  For 
graphical  clarity,  the  vertical  structures  of  the  tidal  currents 
for  sites  KB1  and  VR1  2  are  not  displayed  together  with 
the  other  sites  because  of  their  different  character  (i.e., 
eccentricity,  strength,  and  orientation).  A  gradual  change  in 
tidal  characteristics  is  resolved  in  the  ADCP  measurements 
from  depth  cells  about  20  m  above  the  bottom  to  the  deepest 
measured  depth  ccl  1 .  For  M 2  (F igurc  2),  the  semi  major  axes 


decrease  toward  the  bottom,  indicating  weakening  of  the 
tidal  current  speeds,  the  semi  minor  axes  increase  toward 
the  bottom,  indicating  broadening  of  the  tidal  current  ellipses 
and  counterclockwise  rotation  of  the  current  vectors  around 
these  ellipses,  and  the  phases  decrease  toward  the  bottom, 
indicating  upper  currents  lagging  bottom  currents.  Changes 
in  the  ellipse  orientation  arc  small,  with  some  clockwise 
rotation  of  the  ellipse  toward  the  bottom.  This  structure 
agrees  with  predictions  from  tidal -forced,  bottom  Ekman 
theory  [e.g.,  sec  Soulsby ,  1983],  as  will  be  further  detailed 
in  section  4  of  this  paper.  The  vertical  structures  of  S2 
currents  (not  shown)  have  the  same  character  as  those  of  M2. 

[15]  The  K\  tidal  currents  (Figure  3)  also  have  gradually 
decreasing  semi  major  axes  and  gradually  increasing 
semi  minor  axes  toward  the  bottom.  However,  although 
the  ellipse  shape  changes  arc  similar  to  those  for  the 
semidiurnal  tides,  the  ellipse  orientation  and  phase  changes 
differ.  The  Kj  tidal  ellipses  rotate  counterclockwise  toward 
the  bottom,  with  larger  rotations  than  the  clockwise  rota¬ 
tions  of  the  semidiurnal  constituents.  The  Kj  tidal  current 
phase  changes  with  depth  are  much  weaker  than  the  semidi¬ 
urnal  phase  changes.  At  most  sites,  the  K\  phases  slightly 
increase  toward  the  bottom  (bottom  currents  lag  upper 
currents).  These  differences  in  character  between  diurnal 
and  semidiurnal  tidal  currents  near  the  bottom  were  theoret¬ 
ically  demonstrated  by  Kundu  et  al.  [1981]  and  are  due  to  the 
fact  that  the  semidiurnal  and  diumal  frequencies  are  respec¬ 
tively  faster  and  slower  than  the  inertial  frequency  at 
latitudes  greater  than  ~30°. 
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semi-major  axis  (cm/s)  semi-minor  axis  (cm/s)  ellipse  orientation  (deg.) 

Figure  3.  As  in  Figure  2  but  for  Kj  tidal  currents. 


[i6]  Figure  4  shows  SS6  tidal  current  ellipses  as  repre¬ 
sentatives  of  typical  tidal  current  structure.  At  middepths, 
the  diurnal  and  semidiurnal  current  ellipses  align  in  direc¬ 
tion  and  remain  relatively  unchanged  over  a  large  depth 
range.  Approaching  the  bottom,  all  the  ellipses  shorten  and 
broaden,  the  diurnal  ellipses  veer  drastically  counterclock¬ 
wise  together,  and  the  semidiurnal  ellipses  veer  slightly 
clockwise.  Phase  changes  are  difficult  to  see  in  this  type  of 
graphic,  but  semidiurnal  currents  strongly  lead  near  the 
bottom,  while  diurnal  currents  weakly  lag. 

[n]  Some  departures  from  these  general  trends  should  be 
noted.  The  Kj  tidal  structure  is  somewhat  different  at  sites 
(not  shown)  K.B1,  VR1,  and  VR2.  At  site  SS2,  there  is  a  Ki 
phase  drop  of  14°  from  15  to  23  m  above  the  bottom.  This 
vertical  anomaly  could  possibly  be  caused  by  the  Po  River 
plume,  following  the  hypothesis  of  Chavanne  et  ai  [2007] 
that  the  plume  caused  an  abrupt  phase  drop  in  surface 
tidal  currents  in  a  horizontal  band  all  along  the  Italian  coast 
as  measured  by  their  high-frequency  radars.  Site  VR4 
shows  an  increase  in  the  phase  toward  the  surface. 

3.  Models 
3.1.  NCOM 

[ i k]  The  Navy  Coastal  Ocean  Model  (NCOM)  is  a  3-D, 
finite-difference  numerical  model  based  on  the  primitive 
equations  and  the  hydrostatic,  Boussinesq,  and  incompress¬ 
ible  assumptions.  NCOM  specifies  vertical  mixing  in  terms 
of  vertical  eddy  coefficients  that  are  calculated  according 
to  either  the  Mellor-Yamada  2.0  or  2.5  turbulence  closure 
scheme.  The  NCOM  runs  for  this  study  used  the  2.0  scheme. 


which  is  described  by  Mellor  and  Yamada  [1974].  Imple¬ 
mentation  of  this  scheme  in  NCOM  and  details  of  the  model 
are  described  by  Martin  [2000]. 

[19]  NCOM  was  implemented  for  the  entire  Adriatic  Sea 
on  a  1.02-km  horizontal  grid.  A  difference  from  previous 
descriptions  of  NCOM  is  that  the  version  of  NCOM  used 
here  was  modified  to  allow  the  use  of  generalized  sigma 
coordinates,  where  the  fractional  sigma  layer  thickness  can 
vary  horizontally  as  well  as  vertically  and  sigma  layers  can 
be  masked  to  land  as  the  bottom  shallows.  The  vertical  grid 
used  for  this  study  was  set  up  with  the  following  properties: 
40  layers  in  deep  water,  with  a  gradual  reduction  to  7  layers 
in  shallow  water,  logarithmic  expansion  of  the  grid  away 
from  the  surface  and  bottom  to  provide  increased  resolution 
in  the  surface  and  bottom  boundary  layers  (with  expansion 
factors  of  1.14  and  1.25,  respectively),  fairly  consistent 
resolution  near  the  surface  everywhere  and  near  the  bottom 
in  water  shallower  than  about  100  m  (with  layer  thicknesses 
at  the  surface  and  bottom  of  about  1.0  and  0.26  m, 
respectively),  nearly  horizontal  layers  in  the  upper  half  of 
the  water  column,  and  less  slope  of  the  layers  than  with  a 
regular  sigma  coordinate  grid  in  the  lower  half  of  the  water 
column  (except  near  the  bottom). 

[20]  The  NCOM  vertical  grid  spacing  for  the  bottom  24  m 
of  the  water  column  was  very  similar  at  all  1 2  JRP  mooring 
locations  where  the  bottom  depths  were  greater  than  25  m 
There  were  between  16  and  18  layers  in  the  bottom  24  m, 
with  vertical  grid  spacing  in  this  region  ranging  from  0  2  m 
to  3.4  m.  The  distance  between  the  bottom  and  the  center  of 
the  deepest  sigma  layer  varied  by  less  than  2  em  from  site  to 
site  and  averaged  0.12  m 
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22  m  3  m 


Figure  4.  Measured  tidal  current  ellipses  (cm/s)  at  station  SS6  for  four  different  depths.  M2  ellipses  are 
drawn  in  magenta,  Ki  ellipses  arc  drawn  in  red,  S2  ellipses  are  drawn  in  yellow,  Oi  ellipses  are  drawn  in 
blue,  Pi  ellipses  arc  drawn  in  green,  N2  ellipses  arc  drawn  in  cyan,  and  K2  ellipses  arc  drawn  in  black. 
Dots  indicate  relative  phasing.  The  directions  of  rotation  of  the  ellipses  with  time  arc  from  the  dots  around 
the  ellipses  to  the  gaps.  Distance  from  the  bottom  to  the  midpoint  of  the  ADCP  depth  cells  arc  60,  41,  22, 
and  3  m  as  labeled. 


[21]  Bottom  stress,  fj,  was  calculated  from  bottom  ve¬ 
locity,  Ub,  using  the  quadratic  law, 

t$  =  pQ|ri;|aV  (1) 

A  logarithmic  velocity  profile  in  the  bottom  layer  was 
assumed  in  NCOM,  and,  thus,  Cj  was  adjusted  for  slight 
spatial  variations  of  the  bottom  layer  thickness  according  to 
the  equation, 


where  k  =  0.4  is  von  Karman’s  constant  and  A zb  is  the 
bottom  layer  thickness.  The  bottom  roughness  length,  z0, 
was  kept  constant  in  NCOM  at  a  value  of  0  003  m.  The 
implementation  of  equation  (2)  produced  values  of  C^at  the 
simulated  JRP  mooring  locations  from  0.01 1  to  0.012. 

[22]  NCOM  was  previously  run  for  the  Adriatic  [Martin 
et  aL,  2006]  to  investigate  the  total  circulation  of  the  Sea 
with  as  realistic  a  simulation  as  possible.  Therefore,  accurate 
simulation  of  the  tides  was  an  important  component  of  this 
goal  as  the  fractional  variance  of  measured  northern  Adriatic 
currents  that  could  be  explained  by  tides  ranged  from  1 0%  to 


64%  [see  Martin  et  aiy  2006,  Table  2],  Tidal  forcing  was 
from  tidal  sea  surface  height  and  depth-averaged  velocities 
that  were  prescribed  at  the  open  boundary  in  the  northern 
Ionian  Sea.  These  values  were  taken  from  the  Oregon  State 
University  Mediterranean  (Oj,  Kj,  M2,  and  S2)  and  global 
(Qi,  Pi,  N2,  and  K2)  tidal  databases.  Tidal  potential  forcing 
was  used  in  the  interior  of  the  model  for  these  eight 
constituents,  with  sensitivity  studies  showing  that  about 
12%  of  the  M2  and  7%  of  the  Ki  tidal  elevations  in  the 
northern  Adriatic  eould  be  explained  by  this  direct  astro¬ 
nomical  forcing  [Martin  et  a! .,  2006].  The  Adriatic  NCOM 
version  used  in  this  paper  included  the  total  circulation  as 
well  as  the  tides,  just  as  in  the  previous  study. 

[23]  NCOM  tidal  solutions  were  extracted  by  harmonic 
tidal  analysis  (Q,,  Oh  P],  N2,  M2,  S2,  and  K2)  using 

model  results  from  l  September  2002  to  29  April  2003. 
Quantitative  evaluation  of  NCOM-simulated  tidal  sea  sur¬ 
face  height  for  the  Adriatic  was  previously  done  by  com¬ 
parison  with  data  from  27  International  Hydrographic 
Organization  stations  [see  Martin  et  al.y  2006,  Tabic  1] 
The  RMS  error  was  lower  than  1.6  cm  in  amplitude  for  all 
constituents  and  was  9  and  7  degrees  in  phase  for  M2  and 
Kj,  respectively.  Although  this  previous  study  was  done 
with  an  Adriatic  version  of  NCOM  with  much  less  vertical 
resolution  near  the  bottom  than  the  simulation  used  here,  the 
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overall  tidal  accuracy  outside  the  bottom  zone  is  expected  to 
be  similar.  Martin  et  ai  [2006]  also  showed  good  qualita¬ 
tive  agreement  between  the  JRP  mooring  tidal  currents  and 
the  NCOM  tidal  currents  but  no  quantitative  comparison 
was  presented 

3.2.  Quoddy 

[24]  “Quoddy”  is  a  finite-element  numerical  model  based 
on  the  3-D  nonlinear  shallow  water  equations  and  using  the 
hydrostatic,  Boussinesq,  and  incompressible  assumptions. 
Quoddy  specifies  vertical  mixing  using  a  Mellor-Yamada 

2.5  turbulence  closure  scheme  [Mellor  and  Yamada ,  1082] 
with  the  improvements  of  Galperin  et  ai  [1988].  Imple¬ 
mentation  of  this  scheme  inside  Quoddy  and  details  of  the 
model  arc  fully  described  by  Lynch  et  ai  [1996]. 

[25]  Quoddy  was  implemented  for  simulation  of  the 
Adriatic  tides  using  a  finite-element  mesh  with  typical  node 
spacing  of  500  m  in  coastal  areas  and  44  km  in  deep  areas. 
Quoddy  used  sigma  vertical  coordinates  with  21  vertical 
layers.  Quoddy  used  a  sinusoidal  spacing  of  vertical  levels, 
with  highest  resolution  in  both  the  surface  and  bottom  layers 
and  lower  resolution  at  middepth.  At  the  JRP  locations, 
between  9  and  1 5  of  these  levels  were  located  in  the  bottom 
24  m,  and  vertical  node  spacings  were  between  1 .0  m  and 

4.6  m  there  Quoddy  uses  the  quadratic  law  for  calculation 
of  bottom  stress,  but  with  a  constant  bottom  drag  coefficient 
Cd  -  0.003  because  the  bottommost  node  is  always  1  m  off 
the  bottom  everywhere.  Evaluating  equation  (2)  with  A zh  = 
2  m  gives  an  NCOM  equivalent  Cj  at  1  m  of  0.0047.  The 
Quoddy  value  for  Cj  was  found  through  a  series  of 
numerical  tuning  experiments.  A  complete  description  of 
the  Adriatic  Quoddy  setup  is  presented  by  Janekovic  and 
Kuzmic  [2005]. 

[26]  Unlike  NCOM,  Quoddy  was  not  forced  from  tidal 
database  values  at  the  boundary.  Instead,  Quoddy  was 
iteratively  coupled  to  a  linear,  3-D,  finite-element  model 
and  inverse  system,  “Truxton/Fundy”  [Lynch  and  Naimie , 
1993],  to  determine  sea  level  boundary  conditions  in  the 
Strait  of  Otranto.  The  Truxton/Fundy  data  assimilation 
system  used  data  from  six  coastal  tide  gauge  stations  to 
produce  optimized  boundary  conditions.  These  boundary 
conditions  were  then  used  in  a  Quoddy  run  to  produce 
residual  errors  at  the  stations,  which  were  then  in  turn  used  in 
Truxton/Fundy  as  data  to  produce  an  update  to  the  boundary 
conditions.  The  procedure  was  iterated  to  obtain  boundary 
conditions  specifically  optimized  for  Quoddy.  Truxton/ 
Fundy  used  A  =  1  *10  5  m/s  as  the  linear  frictional  parameter 
and  a  constant  vertical  viscosity  of  0,04  m2/s.  Boundary 
conditions  for  the  Oj,  Pi,  K|,  N2,  M2,  S2,  and  K2  tides  were 
each  solved  for  separately  and  then  Quoddy  was  run  sepa¬ 
rately  for  each  constituent  and  once  with  all  seven  constitu¬ 
ents  together.  Janekovic  et  ai  [2003]  conducted  numerical 
experiments  with  Quoddy  in  the  Adriatic  using  direct  astro¬ 
nomical  forcing.  They  found  much  less  effect  than  Martin 
et  ai  [2006]  with  less  than  1%  contribution  to  M2  ele¬ 
vation  amplitude  and  6%  contribution  to  K]  elevation 
amplitude  at  the  northwest  end  of  the  Adriatic  Therefore, 
direct  astronomical  forcing  was  not  used  in  later  Quoddy 
runs.  Janekovic  and  Kuzmic  [2005]  provides  further  details 
on  the  data  assimilation  and  tide  forcing  procedures. 

[27]  Quoddy  tidal  solutions  were  extracted  by  harmonic 
tidal  analysis  (Oi,  Pi,  Kj,  N2,  M2,  S2,  and  K2)  using  the 


seven-constituent  model  results  from  1  February  1982  to 
16  May  1982,  with  inference  techniques  to  resolve  the  K^Pj 
and  S2/K2  constituent  pairs.  Janekovic  and  Kuzmic  [2005] 
previously  validated  Quoddy  tidal  simulations  using  coastal 
tide  gauge  station  data,  rotary  current  meter  data  at  eight 
sites,  and  an  ADCP  record.  In  comparisons  with  data  from 
3 1  coastal  tide  gauge  stations,  the  RMS  error  was  lower  than 
0.8  cm  in  amplitude  for  all  constituents  and  was  10  and 

6  degrees,  respectively,  in  phase  for  M2  and  K!  [see 
Janekovic  and  Kuzmic ,  2005,  Table  3].  Comparisons  with 
tidal  velocity  measurements  also  showed  good  agreement, 
especially  with  respect  to  simulation  of  tidal  ellipse  orien¬ 
tations.  The  simulated  tidal  current  variability  with  depth 
qualitatively  matched  the  major  characteristics  observed  in 
ADCP  data  taken  at  the  single  location. 

4.  Theory 

[28]  The  governing  momentum  equation  for  tidal  flow  in 
an  unstratified  environment  over  a  localized  area  is 

du  7  _  -  Of.  du\ 

where  u  is  the  horizontal  velocity  vector,  /  is  time,  f  is  the 
Coriolis  parameter  (/)  times  the  vertical  unit  vector,  g  is 
gravitational  acceleration,  rj  is  sea  surface  elevation,  z  is 
vertical  height  above  the  bottom,  and  Az  is  a  coefficient  of 
eddy  viscosity.  Here  we  have  made  the  hydrostatic 
assumption,  neglected  horizontal  advective  fluxes  and 
horizontal  diffusion  of  momentum,  and  assumed  that  direct 
astronomical  forcing  of  the  tides  is  negligible  compared  to 
the  co-oscillating  tide  because  of  the  localization  of  the  area 
considered.  Because  of  the  strength  of  the  Adriatic  tidal 
currents  and  their  horizontal  spatial  scales,  the  horizontal 
advective  fluxes  are  estimated  to  be  <1%  of  and  the 
horizontal  diffusion  of  momentum  fluxes  are  estimated  at 
0.1%  of  ^  for  the  northern  Adriatic  from  late 

September  to  early  May,  stratification  is  generally  weak 
[sec  Jeffries  and  Lee ,  2007,  Figure  2]  and  a  conservative 
estimate  of  the  Burger  number  over  the  entire  water  column 
for  our  application  is  0.01.  Therefore,  especially  in  the 
bottom  Ekman  layer,  stratification  should  play  a  very  minor 
role  in  the  momentum  balance. 

[29]  The  complication  that  remains  is  specifying  the 
vertical  structure  of  the  eddy  viscosity.  The  easiest  choice 
is  to  assume  that  it  is  constant  through  a  bottom  boundary 
layer.  With  this  choice,  solutions  for  tidal  problems  with 
friction  were  found  by  Sverdrup  [1927],  Munk  et  ai  [1970], 
Kundu  et  ai  [1981  ],  Prandle  [1982],  and  others.  If  the  tidal 
current  above  the  boundary  layer  (where  friction  is  unim¬ 
portant)  for  a  particular  constituent  is  u  =  exp(— iuit)i  + 

7  cxp(— /urt)j},  then  the  tidal  solution  in  the  boundary  layer 
is  u(z)  =  S{i/(z)exp(-zutf)i  +  v'(z)exp(-zu;Oj}  with 


/,  .  («  +  h) 

u  (z)  =  a  exp 

(-1  ±1), 

\f  -  w| 

1  2A:  ' 

(-a  + 17) 

+  - - 2  ~  exp 
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and 


v'(z)  =  7  +  -  exp 


2 

(k*  +  7) 


exp 


(-1  ±0 
(-1  +*) 


I /-<•* 


2^. 


+  uJ 


2A 


(5) 


where  a  and  7  are  complex  constants,  u;  is  the  tidal 
constituent  angular  frequency,  the  positive  signs  in  the  first 
exponential  terms  are  for  the  case  u ;  >/ (semidiurnal  tides  of 
the  Adriatic),  and  the  negative  signs  arc  for  the  case  u ;  </ 
(diurnal  tides  of  the  Adriatic).  The  boundary  conditions  that 
were  used  are  a  no-slip  condition  at  the  bottom  and  a  no-stress 
condition  at  00.  Equations  (4)  and  (5)  are  only  valid  for  the 
Northern  Hemisphere.  Further  sign  changes  are  needed  to 
transform  them  for  the  Southern  Hemisphere. 

[30]  A  best  fit  to  these  equations  was  sought  for  the 

ADCP  data.  Az  was  treated  as  an  unknown  and  a  and  7 
were  determined  from  tidal  velocity  means  from  ADCP 
depth  cells  more  than  24  m  above  the  bottom,  excluding  the 
depth  cell  nearest  the  surface.  Equivalent  Ekman  depths  for 
equations  (4)  and  (5)  differ  for  the  clockwise  (CW)  and 
counterclockwise  (CCW)  rotary  components  of  the  tidal 
flows  (see  section  6),  with  the  CW  rotary  depth,  HE  = 
7Ty/lA2/(\f  —  u>|),  always  greater  than  both  the  steady  flow 
and  CCW  rotary  Ekman  depths.  Using  Az  =  8  •  10  4  m2/s 
(see  section  5.1),  the  maximum  Ht  for  the  JRP  locations 
was  24  m  and  occurred  at  site  SS2  for  the  CW  rotary 
component  of  the  K\  tidal  flows.  This  result,  together  with 
the  observed  structure  of  the  tides  in  Figures  2  and  3, 
suggests  that  frictional  effects  do  not  extend  significantly 
beyond  24  m  above  the  bottom.  The  near-surfaee  depth  cell 
was  also  excluded  as  the  tide  solution  was  sensitive  to 
interpolation  methodology  (see  section  2.1).  Tides  for  sites 
SS2,  VR1,  and  VR2  were  not  fit  to  the  infinite-depth 
equations  (4)  and  (5)  because  the  Ekman  depth  extended 
over  the  entire  measured  water  column  at  these  shallow 
sites.  Az  was  found  for  each  site  and  for  the  M2,  S2,  and  Ky 
tides  separately  by  using  a  range  of  possible  values  from  5 
10  5  nr/s  to  2  10' 3  m2/s  in  steps  of  1  •  10  5  m2/s  and 

computing  a  cost  function  as  the  sum  of  the  squared  errors 
between  the  data  and  the  results  of  equations  (4)  (5)  for  all 
depths  with  ADCP  measurements. 

[31]  A  more  complex  expression  for  Az  is  to  use  a  linear 
approximation  and  set  Az  =  0z  +  k ,  where  0  is  a  constant 
coefficient  with  units  of  velocity  and  k  is  a  constant 
kinematic  viscosity.  Prandle  [1982]  solved  equation  (3) 
applied  to  tides  for  this  case.  Following  Prandle  [1982] 
and  forms  given  by  Boas  [1983],  versions  of  equations  (4) 
and  (5)  were  derived  for  the  ease  of  linear  A7  and  boundary 
conditions  of  no  slip  at  the  bottom  and  no  stress  at  00.  They 
are 


u’(z)  =  «  +  J  [ker(v/|/- w|G(r))  ± «  kei(v/|/  -  w|G(z))] 

+  y  [ker(\/(/  +  u>)G(z))  -  i  kei(\/(/  +  w)G(z))] 

(6) 


and 

v'(z)  =  7  +  Y  [-'  ker(v/|/-u;|G(z))  ±  kci(v/|/  ~  wIG(z))] 

+  Y  [»  ker(\/ (/  +  w)G(z))  +  kei(\/(/  +  u>)G(z))], 

(V 

with 

(-Q  -  i'7)kcr(v/|/-u>|//)  ±  (-7  +  ;V»)kci(  J\f  -  wftf) 
ker2(v/|7^|//)  +kci2(v/|7rr^f//) 

(8) 

_  (~a  +  /7)ker(v/(/  +  u))//)  -  (7  +  /a)kei(v/(/  +  J)H) 

ker2  ( sj{f  4-  uj)H)  +  kei2  ( y/(f  +  u 7)11 ) 

(9) 

-fa  "»> 

0(2)  -  ^jf  +  IP.  (11) 

where  ker  and  kei  are  Kelvin  functions  of  order  zero.  For 
terms  in  equations  (6)  (8)  with  alternative  signs,  the  positive 
signs  arc  for  u  <f  and  the  negative  signs  are  for  u;  > /.  This 
set  of  equations  is  only  valid  for  the  Northern  Hemisphere. 

[32]  As  was  done  for  equations  (4)  and  (5),  a  best  fit 
was  sought  between  the  ADCP  data  and  equations  (6)~(1 1 ). 
An  identical  procedure  to  the  one  described  for  those 
equations  was  used  for  these,  with  0  varied  from  4  10  s  m/s 
to  1.6  •  10  1  m/s  in  steps  of  4  ■  10  5  m/s  and  k  varied 
from  1  •  10  5  m2/s  to  1  •  10  1  m2/s  in  steps  of  1  •  10  s  m2/s. 
0  and  k  were  optimized  jointly. 

[33]  Equations  (6)  and  (7),  using  best  fit  values  for  0  and 
k ,  do  not  converge  asymptotically  to  a  and  7  for  increasing 
z  as  quickly  as  equations  (4)  and  (5)  using  best  fit  values  for 
Az.  This  means  that  the  influence  of  bottom  friction  can  be 
responsible  for  slight  curvatures  in  tidal  structure  far  from 
the  boundary  and  an  infinite-depth  approximation  is  less 
accurate  for  linear  Az  models  than  for  constant  Az  models. 
Practically,  the  tidal  ellipse  parameters  calculated  using 
linear  Az  fits  have  broad  curvatures  above  the  Ekman  layer 
such  that  their  values  slightly  differ  from  “00”  values 
throughout  the  entire  water  column  at  the  JRP  sites.  Wc 
used  an  iterative  approach  to  cope  with  this  problem  and  to 
reduce  the  sensitivity  of  the  0  and  k  fits  to  the  choices  for 
a  and  7.  Once  the  0  and  k  best  fit  values  were  found,  a 
and  7  were  recalculated  so  that  the  solutions  of  equations  (6) 
and  (7)  would  pass  through  the  means  from  the  ADCP  depth 
cells  more  than  24  m  off  the  bottom  (again  excluding  the 
depth  cell  nearest  the  surface)  at  the  midpoint  of  these  depth 
cells  instead  of  at  00.  Then,  new  best  fit  values  for  0  and  k 
were  found  using  these  new  a  and  7  values.  For  consistency 
of  comparison,  an  analogous  iterative  procedure  was  also 
used  in  constant  Az  fits,  even  though  iteration  had  less  impact 
on  these  solutions. 
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Figure  5.  Best  fit  values  and  RMS  errors  for  theoretical  solutions,  (a)  Best  fit,  depth-constant  Az  values 
with  cireles  for  Kb  dots  for  M2,  and  asterisks  for  S2.  Best  fit  (b)  k  and  (e)  ft  values  for  linear  Az  theory  with 
triangles  for  K|,  diamonds  for  M2,  and  stars  for  S2.  (d)  The  square  root  of  the  mean  of  the  sum  of  all  the 
squared  differences  in  u  and  v  between  theory  and  observations  over  the  bottom  24  m  and  over  a  tidal 
period.  The  thin  line  is  for  depth-constant  Az  theory,  and  the  thick  line  is  for  linear  Az  theory  with  the  same 
symbols  as  in  Figures  5a  5c  for  each. 


[34]  Bottom  stress  was  calculated  numerically  for  applica¬ 
tions  of  equations  (6)  (11)  using  a  forward  finite  difference 
approximation  for  the  vertical  derivative  of  velocity  evalu¬ 
ated  at  the  bottom, 


=  pAz(0)—\ 

02  \  2=0 

%  £*{</(  Az)exp(^)i  +  v/(Az)exp(-/utf)j  J.  (12) 

Numerical  tests  showed  that  u  and  v  as  defined  by 
equations  (6)  and  (7)  are  approximately  linear  near  the 
bottom  for  the  parameter  values  we  are  using,  and  Az  is  also 
nearly  equal  to  k.  Therefore,  the  stress  is  nearly  constant 
(as  in  a  viscous  sublayer),  so  evaluation  of  equation  (12) 
within  this  zone  is  insensitive  to  the  value  of  Az.  For  the 
calculations  of  theoretical  bottom  stress  in  this  paper,  a  Az 
of  1  em  was  used. 

5.  Comparison  of  Theory,  Data,  and  Models 

5.1.  Fitting  Bottom  Ekman  Layer  Theory  to  the  Data 

[35]  Figure  5a  shows  the  values  for  Az  in  equations  (4)  and 
(5)  that  gave  the  best  least  squares  fit  to  the  observations.  Fits 


for  the  Ki  tides  at  site  KB1  were  discarded  here  and  in  all 
subsequent  fits  to  theory  because  of  decreasing  errors  toward 
zero  friction  parameter  values.  That  is,  the  best  fit  procedures 
could  not  mateh  the  observed  Ki  vertical  structures  at  site 
KB1,  likely  because  of  its  low  signal  level  and  signal-to- 
noise  ratio  (the  mid-water-eolumn,  semi  major  axis  was  less 
than  2  cm/s)  and  the  best  fits  tended  toward  vertically  uniform 
tides.  The  results  did  not  vary  much  from  constituent  to 
constituent,  which  demonstrates  the  skill  of  equations  (4) 
and  (5)  in  simulating  the  distinctly  different  vertical 
structures  of  the  semidiurnal  and  diurnal  tides.  The  aver¬ 
age  Az  values  were  7.1  •  10  4  m2/s,  7.8  •  10  4  m2/s,  and 
8.4  10  4  m2/s  for  K|,  M2,  and  S2,  respectively,  and  7.8  • 

10  4  m2/s  overall.  Variation  from  site  to  site  was  also  weak, 
with  marginally  higher  values  at  sites  SS4  5. 

[36]  The  thin  lines  in  Figure  5d  show  the  residual  errors 
between  the  best  fit  constant  Az  theory  and  the  data  from  the 
bottom  24  m  of  the  water  column.  The  average  error  was 
0  31  cm/s,  with  higher  average  error  (0.40  em/s)  for  the 
stronger  M2  tides  compared  to  the  S2  (0.25  em/s)  and  Kj 
(0.29  em/s)  tides.  Overall,  the  errors  peaked  at  0  61  cm/s  for 
the  M2  tides  at  site  SS10.  The  and  S2  best  fit  solutions 
both  had  error  peaks  at  site  SS4. 

[37]  By  factoring  out  /?,  the  linear  form  for  Az  becomes 
/ 3(z  +  k/fl),  where  kJ/3  is  a  scale  height.  This  term  eannot 
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simply  be  neglected  because  kcr  and  kei  approach  infinity 
as  z  approaches  zero,  but  Prandle  [1982]  assumes  that  it  is 
vanishingly  small.  This  is  in  accord  with  the  classic 
“logarithmic  layer”  approach  of  assigning  Az  =  ku+z  +  v 
[Wimbush  and  Munk ,  1970],  where  w*  is  the  friction 
velocity  and  is  equal  to  the  square  root  of  the  bed  shear  stress 
divided  by  density  and  1/  =  1 .2  •  1 0  6  m2/s  is  the  molecular 
kinematic  viscosity.  Soulsby  [1983]  uses  this  form,  addition¬ 
ally  neglecting  the  scale  height  term,  now  iV(kw*)>  by  setting 
the  condition  of  no  slip  at  z  =  Zo  ( Soulsby  [1983]  used  z0  = 
0.0009  m)  instead  of  z  =  0.  Therefore,  as  a  first  attempt  to 
improve  upon  the  results  of  optimizing  equations  (4)  and  (5), 
equations  (6)  (11)  were  used,  k  was  set  equal  to  vy  and  a 
best  fit  to  the  data  was  sought  by  optimizing  0  -  ku *. 

[38]  This  form  for  Az  was  expected  to  be  a  more  realistic 
approximation  than  a  constant  Az.  However,  when  using 
molecular  kinematic  viscosity  and  solving  for  an  optimal  w* 
on  the  basis  of  the  ADCP  data,  this  did  not  prove  to  be  true. 
The  RMS  error  in  the  bottom  24  m  remaining  after  optimi¬ 
zation  (not  shown)  was  lower  in  only  11%  of  the  fits 
compared  to  using  constant  Az.  Despite  being  worse  than 
constant  Az  fits  in  a  large  majority  of  cases,  th cAz  =  ku*z  +  v 
fits  were  better  in  simulating  the  measured  semi  major  axis 
amplitudes  in  74%  of  the  cases.  The  reason  that  this  result 
did  not  translate  to  overall  better  performance  was  that  the 
Az  =  ku+z  +  v  form  better  simulated  semi  -minor  axis 
amplitudes  in  only  46%  of  the  cases,  ellipse  orientation  in 
only  6%  of  the  cases,  and  phase  in  only  1 1%  of  the  cases. 

[39]  Because  of  this  result,  the  fits  were  redone  allowing 

both  0  and  k  to  vary  as  specified  in  section  4.  For  the  S2  fit 
at  site  SS4,  k  was  varied  to  higher  values  (2  10  1  m2/s)  in 

order  to  find  its  best  fit  value  close  to  the  maximum  k 
considered  for  other  sites.  Figures  5b  and  5c  show  the 
resulting  optimum  values  for  these  parameters.  Values  for 
Ki  at  sites  VR4  5,  for  S2  at  site  SS6,  and  for  all  constit¬ 
uents  at  site  CP2  are  not  shown  because  the  fit  selected  the 
lowest  value  of  0  (4  ■  10  5  m/s)  and  thus  a  true  minimum 
was  not  found.  Clearly  for  these  cases,  the  best  fit  was 
approaching  the  constant  Az  case,  which  is  mathematically 
equivalent  to  using  equations  (6)  (1 1 )  with  0  =  0.  For  other 
cases,  the  optimum  0  occupied  a  fairly  small  range  of  values 
independent  of  coefficient  and  site.  Values  for  Ki  and  S2  at 
site  SS4,  for  Kj  at  site  SS9,  and  for  M2  at  site  SS10  were 
somewhat  larger,  with  values  more  than  two  times  the 
overall  average.  Average  0  values  were  3.2  •  10  4  m/s, 
2.7-10  4  m/s,  and  2.7  •  10  4  m/s  for  K,.  M2,  and  S2 
respectively  and  2.9  •  10  4  m/s  overall. 

[40]  Optimum  k  values  (Figure  5b)  were  more  than  100 

to  more  than  800  times  greater  than  molecular  kinematic 
viscosity.  Optimum  k  also  occupied  a  fairly  small  range  of 
values  independent  of  coefficient  and  site,  with  a  weak 
trend  toward  lower  values  from  southwest  to  northeast 
along  the  SS  line.  Average  k  values  were  4.6  10  4  m2/s, 

4.7  10  4  m2/s,  and  5.7  10  4  m2/s  for  K,,  M2,  and  S2 

respectively  and  5,0  •  10' 4  m2/s  overall. 

[41]  The  thick  lines  in  Figure  5d  show  the  residual  errors 
between  the  best  fit  linear  Az  theory  and  data.  We  have 
included  the  errors  for  fits  with  the  minimum  0  for 
comparison  to  the  constant  Az  residual  errors.  Linear  Az 
theory  has  lower  errors  than  constant  Az  theory  in  91%  of 
the  cases.  The  average  error  decreased  by  30%  to  0.22  cm/s. 
Average  errors  for  particular  constituents  were  0.22  cm/s  for 


Ki,  0.26  cm/s  for  M2,  and  0.19  cm/s  for  S2.  Comparing  the 
error  for  each  tidal  current  parameter,  linear  Az  better 
simulated  the  semi  major  axis  in  91%  of  the  cases,  the 
semi -minor  axis  in  74%  of  the  cases,  the  ellipse  orientation 
in  51%  of  the  cases,  and  the  phase  in  83%  of  the  cases. 
Considering  only  the  ellipse  orientation  for  K!  and  only  the 
phase  for  M2  and  S2  (i.e.,  the  respective  dominant  angle 
changes),  linear  Az  performed  better  in  83%  of  the  cases 

[42]  Figures  6  and  7are  representative  of  the  fits.  For  both 
diurnal  and  semidiurnal  tides,  constant  Az  theory  (black 
curves)  produces  a  distinct  bulge  in  the  semi  major  axes 
and  semi  minor  axes  at  the  top  of  the  bottom  Ekman  layer 
analogous  to  the  increase  in  current  speed  produced  for  a 
steady,  constant  Azy  bottom  Ekman  layer.  However,  such  a 
bulge  is  generally  not  seen  in  the  ADCP  data  and  a  linear  Az 
(red  curves)  is  better  able  to  match  the  currents  in  this  area 
because  it  also  lacks  a  distinct  bulge  and  instead  has  a 
broad  maximum.  Both  constant  Az  and  linear  Az  theory 
match  well  the  sharp  changes  in  diumal  ellipse  orientation 
and  semidiurnal  phase  observed  in  the  bottom  ADCP  depth 
cells,  but  constant  Az  theory  generally  overpredicts  the 
change  around  the  depth  level  of  the  bottom  depth  cell.  Note 
that  it  is  the  constant  component  of  linear  Az  theory  that  is 
primarily  responsible  for  reproducing  these  angular  changes. 
Although  linear  Az  theory  matches  the  semi  major  axis, 
semi  minor  axis,  diumal  ellipse  orientation,  and  semidiurnal 
phase  changes  more  accurately  than  constant  Az  theory  in 
most  cases,  weak  changes  in  the  semidiurnal  ellipse  orienta¬ 
tion  and  in  diumal  phase  match  constant  Az  theory  better  in 
half  the  cases. 

[43]  Sensitivity  studies  were  performed  to  investigate  the 
impact  on  the  best  fit  values  from  using  different  cost 
function  forms.  Specifically,  20  m  and,  at  some  sites,  30  m 
were  tried  as  selection  heights  for  the  “Ekman  layer”  and 
calculation  of  a  and  7.  Also,  for  these  fits  (and  an  additional 
fit  with  a  24-m  Ekman  layer)  the  cost  functions  were 
evaluated  only  over  the  bottom  layer  and  not  over  the  entire 
water  column,  and  no  iteration  of  a  and  7  was  performed. 
The  results  were  relatively  insensitive  to  these  variations  in 
methods  and  cost  function  definitions.  Average  /?  and  k  best 
fit  values  (not  including  three  outliers)  for  the  20-m  cost 
function  were  3.3  •  10  4  m/s  and  4.2  •  10  4  m2/s,  respec¬ 
tively.  The  average  constant  Az  value  for  this  cost  function 
was  8.0  •  10~4  nr/s. 

5.2.  Tidal  Current  Vertical  Structure  Comparison 

[44]  The  degree  of  mismatch  between  the  observations 
and  model  simulations  in  the  bottom  Ekman  layer  is 
sensitive  both  to  how  the  model  represents  the  bottom  layer 
and  to  the  model  solution  errors  above  the  Ekman  layer. 
Separation  of  these  two  effects  is  desirable  but  cannot  be 
simply  achieved.  The  approach  we  have  taken  is  to  estimate 
biases  in  the  upper  water  column  and  remove  these  from  the 
entire  water  column  before  calculating  error  statistics  in  the 
Ekman  layer.  This  is  an  imperfect  solution  since,  using 
linear  Az  theory  for  an  example,  it  removes  bias  from  the 
a  and  7  terms  in  equations  (6)  and  (7)  but  does  not  remove 
the  effect  of  model  solution  error  in  B>  and  Df  caused  by 
incorrect  a  and  7  values  in  equations  (8)  and  (9).  Thus  the 
depth  structure  of  the  tides  depends  on  a  and  7  in  a  complex 
way.  Sometimes  removal  of  the  upper  water  column  bias  will 
cause  error  in  the  Ekman  layer  to  increase.  However,  all 
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Figure  6.  Observations,  theoretical  fits,  and  model  simulations  of  tidal  currents  at  SS6.  The  four 
panels  show  tide  ellipse  parameters  versus  depth  as  individually  labeled  The  blue  curves  are  measured 
values  with  dashed  lines  indicating  95%  confidence  limits,  the  black  curves  are  from  the  constant  Az  best 
fit,  the  red  curves  are  from  the  linear  Az  best  fit,  the  green  curves  are  from  NCOM,  and  the  magenta 
curves  are  from  Quoddy.  The  model  simulations  have  been  interpolated  to  ADCP-observed  depth  levels 
and  the  upper  water  column  biases  have  been  removed.  Note  that  the  x  axis  scale  for  ellipse  orientation 
has  a  larger  range  than  the  x  axis  scale  for  phase  because  the  dominant  diurnal  angular  response  is  in 
ellipse  orientation. 


errors  were  calculated  with  and  without  bias  removal  and, 
overall,  the  errors  are  lower  when  the  bias  is  removed.  Bias 
removal  increases  the  NCOM  average  (over  sites)  RMS 
error  in  the  Ekman  layer  for  the  Ki  tide,  but  decreases  this 
RMS  error  for  all  the  Quoddy  tides  and  for  the  M2  and  S2 
NCOM  tides.  Therefore,  to  better  examine  how  the  models 
represent  the  Ekman  layers,  errors  with  upper  water  column 
bias  removed  are  the  ones  that  are  used. 

[45]  Figures  8  10  show  the  calculated  upper  water 
column  biases  from  NCOM  and  Quoddy,  These  were 
calculated  as  the  difference  in  ellipse  parameters  using 
average  a  and  7  values  calculated  from  the  models  and  data 
over  the  same  depth  ranges  determined  from  ADCP  depth 
cells  higher  than  24  m  off  the  bottom.  Quoddy  tends  to  have 
semi  major  axes  that  are  biased  too  large  and  removal  of 
these  biases  greatly  helps  to  lower  RMS  errors  for  the  Ekman 
layer.  NCOM  M2  and  S2  tidal  currents  all  lead  the  observed 
tidal  currents  except  at  site  KB1.  Removal  of  this  phase  bias 
greatly  lowers  the  NCOM  Ekman  layer  errors.  The  Ki  biases 
for  the  semi  major  axes  and  semi  minor  axes  are  especially 
interesting  as  both  NCOM  and  Quoddy  have,  uniquely  for 
these  cases,  very  similar  patterns  of  bias.  This  suggests  the 
presence  of  a  dynamics  or  parameter  error  that  is  common  to 


both  models  and  affects  the  Ki  tidal  currents  without  causing 
large  effects  in  the  M2  or  S2  tidal  currents. 

[46]  Figure  11  shows  the  total  RMS  error  for  each  tidal 
constituent  evaluated  for  the  layer  within  24  m  of  the  bottom. 
The  bias  correction  described  in  the  previous  paragraphs 
has  been  applied.  The  error  from  the  best  fit  for  linear  Az 
from  Figure  5d  is  also  shown  for  comparison.  The  average 
errors  for  NCOM  are  0.37  em/s  for  K,,  0.36  cm/s  for  M2, 
and  0.33  cm/s  for  S2  (0.35  cm/s  overall)  and  for  Quoddy  are 
0.33  em/s  for  Kj,  0.44  cm/s  for  M2,  and  0.30  cm/s  for  S2 
(0.36  em/s  overall).  Results  for  linear  Az  theory  averaged 
38%  lower  than  these  values.  NCOM  and  Quoddy  have 
similar  error  levels,  but,  near  1  stria,  where  M2  tidal  currents 
are  especially  strong  and  the  observed  tidal  structure 
changes  strongly  with  depth,  NCOM  matches  the  currents 
exceptionally  well  while  Quoddy  has  only  weak  tidal 
current  depth  changes.  Although  linear  Az  theory,  NCOM, 
and  Quoddy  all  match  the  observed  tidal  structure  in  the 
bottom  Ekman  layer  relatively  well,  better  agreement  could 
be  achieved  before  reaching  the  error  level  of  the  measure¬ 
ments,  as  directly  illustrated  in  Figures  6  and  7  for  site  SS6. 

[47]  Figures  12  14  show  the  breakdown  of  the  Ekman 
layer  RMS  error  according  to  tidal  ellipse  parameter.  NCOM 
averaged  the  lowest  errors  for  M2  semi  minor  ellipse  axes, 
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F  igure  7.  As  in  Figure  6  but  for  M2  tidal  currents  at  SS6.  Note  that  the  x  axis  seale  for  phase  has  a  larger 
range  than  the  x  axis  scale  for  ellipse  orientation  because  the  dominant  semidiurnal  angular  response  is  in 
phase. 


but  linear  Az  theory  averaged  the  lowest  errors  for  all  other 
constituents  and  ellipse  parameters.  Quoddy  simulated  the 
M2  and  S2  semi  major  axis  changes  relatively  well,  with 
only  5%  higher  average  errors  than  linear  Az  theory  but  this 
result  strongly  depends  on  the  removal  of  the  Quoddy 
upper  water  column  bias  (e.g.,  see  Figure  10).  In  contrast, 
the  removal  of  the  upper  water  bias  greatly  increases  the 
Ekman  layer  RMS  errors  for  the  K|  semi  major  axis  at 
sites  SS4  5  for  NCOM  and  site  SS4  for  Quoddy  (i.e,,  the 
three  highest  errors  for  this  parameter  in  Figure  12),  but 
average  (over  stations)  semi  major  axis  Ekman  layer 
RMS  errors  arc  reduced  by  bias  removal.  There  is  a  distinct 
pattern  in  the  degree  that  linear  Az  theory  averaged  lower 
errors  for  the  different  ellipse  parameters,  with  overall  22% 
and  17%  lower  errors  for  semi  major  ellipse  axes  and 
semi  “minor  ellipse  axes,  but  47%  and  53%  lower  errors 
for  ellipse  orientation  and  phase.  In  particular  for  the 
semidiurnal  tides,  linear  Az  theory  averaged  62%  lower 
errors  in  phase. 

5.3.  Time  and  Vertical  Structure  of  Eddy  Viscosity 

[48]  The  analytical  solutions  of  section  4  assume  very 
simple  vertical  structures  for  Az  and  also  assume  that  Az  is 
constant  in  time.  However,  both  NCOM  and  Quoddy  fully 
simulate  Az  with  complex  vertical  and  temporal  structures 
using  Mellor-Yamada  2.0  and  2.5  turbulence  closure 
schemes,  respectively.  These  temporal  and  depth  structures 
were  extracted  from  the  models  for  analysis  and  comparison. 
Az  was  saved  at  hourly  intervals  along  with  other  variables 


for  the  main  run  of  NCOM.  Az  was  derived  using  other  saved 
turbulence  parameters  for  Quoddy  from  a  dedicated  172-day 
run  made  without  nodal  modulation. 

[49]  Figure  15  shows  the  site-averaged  power  spectral 
density  (psd)  of  Az  for  the  bottom  depth  levels  of  NCOM 
and  Quoddy.  They  were  calculated  from  the  Az  time  series 
with  their  means  removed,  using  Welch’s  averaged  periodo- 
gram  method  over  block  lengths  of  1024  h  (~43  days), 
with  50%  overlapping  Hanning  windows  and  no  detrending. 
The  frequency  structures  of  the  NCOM  and  Quoddy  results 
are  similar,  with  psd  peaks  clustered  around  particular 
frequency  bands  near  0,  1,2,  3,  and  4  cycles  per  day  (epd). 
Quoddy  often  has  higher  psd  in  these  bands  and  has  further 
psd  peaks  at  higher  frequencies  that  are  not  excited  in 
NCOM.  In  contrast,  NCOM  has  higher  psd  between  bands 
as  would  be  expected  from  a  model  simulating  the  total 
circulation  rather  than  only  tides.  Both  in  NCOM  and 
Quoddy,  strong  peaks  oceur  at  particular  frequencies 
(marked  with  dotted  lines  in  Figure  15)  within  each  cluster, 
and  each  of  these  peaks  can  be  related  to  specific  frequencies 
of  tidal  constituents  or  tidal  interactions  as  annotated. 

[50]  Figure  16  shows  the  time  means  of  Az  in  the  Ekman 
bottom  layer  for  NCOM  (red)  and  Quoddy  (blue)  at  select 
sites.  NCOM  has  higher  gradients  and  higher  values  of  time 
mean  Az  than  Quoddy  throughout  this  layer.  Both  NCOM 
and  Quoddy  have  considerable  curvature  in  their  time  mean 
Az  depth  structures,  but  the  curvatures  in  the  NCOM  struc¬ 
tures  tend  to  occur  further  off  the  bottom  or  at  higher  Az 
values  and  are  therefore  not  as  evident  in  the  SS8,  CP3,  and 
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Figure  8.  Estimated  model  K!  bias  for  currents  from  heights  more  than  24  m  off  the  bottom.  Squares  are 
for  NCOM,  and  triangles  are  for  Quoddy.  Positive  ellipse  orientation  bias  indicates  that  the  model  tidal 
ellipse  is  rotated  counterclockwise  with  respect  to  the  observed  tidal  ellipse.  Positive  phase  bias  indicates 
that  the  model  tidal  currents  lag  the  observed  tidal  currents.  Mooring  names  arc  shortened  for  ease  of 
display. 


VR6  panels.  The  theoretical  fits  (orange,  black,  and  magenta 
lines)  show  considerable  spread  at  these  scales  from  each 
other  and  from  site  to  site,  but  generally  have  smaller 
gradients  close  to  the  bottom  than  those  of  the  Quoddy  or 
NCOM  Az  time  means. 

[51]  Because  of  the  time  dependence  of  Azy  there  are 
frictional  impacts  at  all  depths  not  described  by  these  time 
means.  Because  the  depth  structure  of  Az  and  the  depth 
structure  of  both  change  with  time  and  with  respect  to 
each  other,  their  interaction  in  equation  (3)  will  vary  in  time 
and  we  cannot  separate  their  frictional  impacts  from  each 
other.  One  way  the  effect  of  Az  can  be  roughly  estimated  is  to 
calculate  a  gain  in  the  semi  major  ellipse  axis  at  tidal 
frequencies  due  to  the  action  of  Az.  This  was  done  by  dividing 
the  semi  major  axis  values  of  harmonically  analyzed  kine¬ 
matic  shear  stress,  (Ajfg),  by  the  semi  major  axis  values  of 
harmonically  analyzed  velocity  gradient,  The  results  for 
Kj  (green)  and  M2  (cyan),  and  for  NCOM  (solid)  and 
Quoddy  (dashed)  are  shown  in  Figure  16.  For  Quoddy, 
these  gains  arc  always  higher  than  the  Quoddy  time  mean 
Az  values  and  often  approach  the  NCOM  time  mean  values 
and  structure  at  depths  close  to  the  bottom.  In  contrast,  the 
NCOM  gains  are  only  slightly  higher  than  the  NCOM  time 
mean  values  for  most  sites  near  the  bottom,  and  they  abruptly 
shift  farther  up  the  water  column  to  lower  and  nearly  depth- 
constant  values.  Above  the  shift,  the  NCOM  gains  are  noisy 
because  of  less  well  determined  semi  major  axis  values  of 


3®.  For  graphical  clarity,  gains  where  the  estimated  error  were 
>50%  of  the  estimated  value  for  cither  axis  arc  not  plotted  in 
Figure  16. 

[52]  Although  the  “impact"  of  Az  alone  can  only  be 
approximated,  the  combined  action  of  Az  and  i.c.  the 
frictional  kinematic  shear  stress,  for  tidal  momentum  at 
specific  frequencies  can  be  exactly  calculated  and  analyzed. 
Here  we  exploit  the  fact  that  the  Fourier  Transform  of  the 
multiplication  of  two  time-dependent  signals  is  the  convo¬ 
lution  integral  of  the  Fourier  Transforms  of  the  individual 
signals.  Therefore,  for  a  particular  tidal  frequency  of  interest, 
the  Fourier  coefficient  of  the  frictional  kinematic  shear 
stress,  ru,  will  be 

Fr\(+u>)  -  f  FAt -  xfidxlu  ( 13) 

J -oc 

where  F  represents  Fourier  transformation. 

[53]  By  this  expression,  if  cither  Az  or  have  dominant 
spectral  peaks,  then  only  a  few  Fourier  coefficients  would 
determine  the  frictional  kinematic  shear  stress  that  impacts 
a  particular  tidal  constituent’s  momentum.  The  convolution 
integral  shifting  property  causes  the  time  mean  of  Az  to 
interact  with  the  ±u>  peaks  of  7^,  and  unless  has  high  psd 
at  other  frequencies,  which  when  shifted  correspond  to 
other  peaks  in  Az ,  then  the  frictional  kinematic  shear  stress 
for  that  tidal  constituent  will  be  mainly  determined  by  the 
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Figure  9.  As  in  Figure  8  but  for  M2  tidal  currents. 


S2  semi-major  axis  S2  semi-minor  axis 


S2  ellipse  orientation  S2  phase 


Figure  10.  As  in  Figure  8  but  for  S2  tidal  currents. 
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Figure  11.  Total  RMS  error  for  all  tidal  currents  from  heights  24  m  or  closer  to  the  bottom.  Currents 
were  corrected  for  upper  water  column  bias  (Figures  8  10)  before  computation  of  RMS.  RMS  errors  for 
(top)  Ki,  (middle)  M2,  and  (bottom)  S2  arc  shown.  Squares  are  from  NCOM,  triangles  are  from  Quoddy, 
and  asterisks  are  from  best  fit  linear  Az  theory. 


time  mean  of  Az.  This  is  the  ease  for  near  bottom  frictional 
kinematic  shear  stresses  in  NCOM  and  Quoddy  and  is  why 
the  gains  track  the  time  means  in  Figure  16. 

[54]  Figure  17  shows  this  in  a  different  form.  Here,  the 
cumulative  sums  of  the  discrete  version  of  equation  (13)  for 
the  direction  of  strongest  are  plotted  for  NCOM  and 
Quoddy  at  M2  and  K\  frequencies  and  at  1.5  and  15  m 
height  above  the  bottom.  Jumps  in  the  cumulative  sums  at 
particular  frequencies  indicate  their  relative  importance  in 
determining  frictional  kinematic  shear  stress.  The  largest 
jumps  are  all  at  zero  frequency.  Because  NCOM  simulates 
the  strong  Western  Adriatic  Current  flow  at  SS4  and  the 
time  mean  of  ^  interacts  in  the  convolution  integral  with 
the  tidal  peaks  of  Az,  there  are  also  large  jumps  at  2  epd 
for  M2  and  1  epd  for  Kj.  Quoddy  completely  lacks  these 
jumps  because  it  does  not  simulate  the  mean  currents.  In  the 
upper  Ekman  layer,  the  Quoddy  kinematic  shear  stresses  are 
weaker  but  the  frequency  contribution  structure  is  similar  to 
that  near  the  bottom.  This  is  not  true  for  NCOM,  where  the 
higher  Az  psd  between  bands  at  lower  frequencies  often 
interacts  with  various  shears  to  interfere  destructively  and 
thereby  reduce  the  frictional  kinematic  shear  stress  that 
would  have  otherwise  been  established  by  a  time  constant  Az. 
This  likely  explains  the  abrupt  shift  and  noise  noted  earlier 
in  the  NCOM  gains  (Figure  16). 

5.4.  Bottom  Stress  Comparison 

[55]  Bottom  stress  at  the  mooring  sites  can  be  calculated 
from  each  of  the  models  and  from  the  best  fit  linear  Az 


theory  results.  Figure  18  shows  examples  of  this  calculation 
for  the  M2  and  K\  tides  at  moorings  SS8  and  VR6.  The 
theoretical  kinematic  bottom  shear  stresses  (solid  black 
curves)  were  calculated  from  equation  (12)  using  the  best 
fit  0  and  k  values,  and  thus  the  time  evolution  of  the  shear 
stresses  takes  the  form  of  an  ellipse  that  is  shaped,  rotated, 
and  phased  in  accord  with  the  near-bottom  tidal  currents. 
The  black  dashed  lines  show  bottom  shear  stresses  that 
would  be  calculated  using  theoretical  values  at  1-m  height 
in  equation  (1)  with  Q  optimized  to  produce  the  same 
maximum  shear  stresses  as  the  solid  black  ellipses.  Use  of  a 
quadratic  drag  law  on  currents  tracing  an  ellipse  in  time 
leads  to  a  double  lobed  bottom  shear  stress  structure;  that  is, 
the  shear  stresses  depart  from  an  elliptical  shape  because  of 
nonlinearity.  Also,  the  orientation  and  phase  of  the  quadratic 
drag  law  stresses  will  match  the  currents  at  the  level  where 
they  are  evaluated  and  miss  any  changes  taking  place  deeper 
in  the  water  column.  The  values  for  Cj  needed  to  match  the 
maximum  shear  stress  magnitudes  of  equation  ( 1 2),  averaged 
over  all  the  theoretically  fitted  sites,  were  the  seemingly  high 
values  of  0.022  and  0.051  for  M2  and  Kj,  respectively. 

[s6]  Also  shown  in  Figure  18  are  the  bottom  shear 
stresses  from  three  models.  Generally,  NCOM  and  Quoddy 
produce  stresses  that  are  similar  to  each  other,  but  smaller 
than  the  theoretical  result.  Both  have  double  lobed  shapes 
from  the  use  of  the  quadratic  drag  law.  Except  for  the  K! 
shear  stresses  at  site  SS4,  NCOM  has  larger  maximum  M2 
and  Ki  bottom  shear  stresses  than  Quoddy  at  all  sites. 
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Figure  12.  K!  tidal  ellipse  parameter  RMS  errors  for  all  currents  from  heights  24  m  or  closer  to  the 
bottom.  Currents  were  corrected  for  upper  water  column  bias  (Figure  8)  before  computation  of  RMS. 
Squares  arc  from  NCOM,  triangles  are  from  Quoddy,  and  asterisks  are  from  best  fit  linear  Az  theory. 

M2  semi-major  axis  M2  semi-minor  axis 


M2  ellipse  orientation  M2  phase 


Figure  13.  As  in  Figure  12  but  for  M2  tidal  currents. 
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S2  semi-major  axis  S2  semi-minor  axis 


Figure  14.  As  in  Figure  12  but  for  S2  tidal  currents.  The  off-scale  RMS  errors  for  KB1  ellipse 
orientation  and  phase  are  16.5°  and  15.4°,  respectively. 


However,  in  19  out  of  these  22  cases,  the  percentage 
increase  is  less  than  what  would  be  expected  simply  because 
of  the  effectively  58%  larger  (1-m)  bottom  drag  coefficient 
in  NCOM.  Thus  this  suggests  that  NCOM  tidal  currents  are 
somewhat  weaker  at  that  level  than  Quoddy  tidal  currents. 
The  percentages  by  which  NCOM  bottom  shear  stresses 
exceed  those  of  Quoddy  generally  increase  toward  the 
northwest  end  of  the  basin,  with  averages  for  the  various 
mooring  lines  of  22%  for  SS,  39%  for  CP,  and  54%  for  VR. 
Note  that  for  both  NCOM  and  Quoddy,  the  actual  bottom 
shear  stresses  associated  with  the  tides  will  be  higher  than 
the  values  in  Figure  18  since  the  double  lobed  ellipses  were 
calculated  using  single  tidal  constituent  velocities  only  and 
the  nonlinearity  of  equation  (1)  will  produce  cross  terms 
with  other  tides  for  Quoddy  and  with  other  tides  and  other 
currents  for  NCOM. 

[57]  The  blue  ellipses  in  Figure  18  are  from  a  depth- 
averaged,  linear,  shallow  water  equation  tide  model  using 
strong-constraint  variational  data  assimilation  of  the  tidal 
observations  in  this  paper.  The  model  is  described  by  Griffin 
and  Thompson  [1996]  and  the  setup  for  the  northern 
Adriatic  is  described  in  Book  [2007].  The  linear  friction 
parameter  was  varied  to  improve  agreement  with  the  data 
and  5  •  10  4  m/s  was  determined  as  the  optimal  value  for 
this  model.  In  the  linear,  depth-averaged  model  equations, 
bottom  stress  is  simply  this  friction  parameter  multiplied  by 
the  depth-averaged  tidal  currents.  Because  of  this,  the 
bottom  stresses  for  a  particular  tidal  constituent  trace  out  a 
true  ellipse.  Also,  since  the  bottom  stresses  are  based  on 


depth-averaged  currents,  they  account  for  very  little  of  the 
broadening,  rotating,  and  phase  shifting  that  takes  place  in 
near-bottom  tidal  currents.  Therefore,  they  arc  much  too 
weak  in  the  cross-axis  direction  of  the  Adriatic.  Contrast¬ 
ingly,  the  optimized  stresses  for  this  model  seem  much  too 
large  in  the  along-axis  direction  compared  to  other  bottom 
shear  stress  estimates. 

6.  Discussion 

[58]  Our  observations  of  near-bottom  tides  show  common 
general  characteristics  relative  to  mid-watcr-column  tidal 
currents:  (1)  tidal  ellipses  shorten  and  broaden,  (2)  diumal 
tidal  current  ellipse  orientations  rotate  strongly  CCW,  and 
(3)  semidiurnal  current  phases  advance  strongly.  All  these 
can  be  qualitatively  explained  by  considering  the  rotary 
components  of  the  tidal  flow  separately  and  examining  each 
of  them  in  a  frame  of  reference  that  is  rotating  such  that  the 
flow  becomes  stationary.  In  these  reference  frames  the  steady 
current  Ekman  solution  above  a  rigid  surface  [ Kundu ,  1990] 
applies,  and  the  speed  and  angular  changes  of  steady  currents 
with  depth  are  equivalent  to  the  amplitude  and  phase  changes 
of  the  respective  rotary  components.  The  steady  current 
Ekman  solutions  are  in  Earth  coordinates  and  are  already  in 
a  rotating  coordinate  system.  Therefore,  additional  rotation 
at  ±uj  to  match  the  tidal  rotary  component  rotation  only  acts 
to  increase  or  decrease  the  equivalent  “Coriolis  parameter.” 

[59]  Thus,  for  the  CCW  rotating  Northern  Hemisphere, 
CCW  rotary  components  have  an  Ekman  depth  of 
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Figure  15.  Average  power  spectral  density  of  Az  for  all  12  JRP  mooring  locations  where  the  bottom 
depths  were  greater  than  25  m.  The  thick  lines  are  calculated  from  the  Quoddy  model  (at  1  m  off  the 
bottom),  and  the  thin  lines  are  calculated  from  the  NCOM  model  (at  ~0.26  m  off  the  bottom).  The  dotted 
lines  mark  specific  frequencies  of  tidal  importance  as  annotated. 


7Ty/2 Az  j(f  +  oj),  CW  rotary  components  have  an  Ekman 
depth  of  7T\/Z4r/(|/  —  u;|),  and  both  rotary  components 
phase  angles  rotate  CCW  with  depth.  However  if  lj  >  f 
(semidiurnal  tides  of  the  Adriatic),  the  reference  frame  itself 
is  rotating  CW  for  CW  rotary  components  and  “Southern 
Hemisphere”  steady  current  Ekman  solutions  must  be  used, 
causing  their  phase  angles  to  rotate  CW  instead  of  CCW. 

[f>o]  For  near-reversing  tides  in  the  Adriatic,  mid-water- 
eolumn  tides  have  rotary  components  that  are  nearly  equal 
in  amplitude  (required  to  produce  reversing  tides)  and  have 
distinct  times  and  orientations  when  the  rotary  veetors  align 
in  direction  producing  maximum  currents  and  setting  the 
ellipse  orientations  and  phases.  The  characteristics  (1  3) 
outlined  at  the  beginning  of  this  section  can  be  explained 
as  follows.  Because  of  their  larger  Ekman  depths,  the  CW 
rotary  component  amplitudes  have  decay  scales  which 
extend  further  up  in  the  water  column  than  the  CCW 
amplitudes.  Consequently  (1)  the  tidal  ellipses  broaden  as 
they  shorten,  and  tidal  current  veetors  rotate  CCW  with 
time  near  the  bottom.  For  diurnal  tides,  both  CW  and  CCW 
rotary  component  have  phase  angles  that  rotate  CCW  with 
depth  (albeit  at  slightly  different  rates)  and,  because  they 
rotate  together,  (2)  changes  in  timing  of  the  maximum 
current  (phase)  are  small  but  the  ellipse  orientations  rotate 
CCW  along  with  the  rotary  component  phases.  For  semi¬ 
diurnal  tides,  CCW  rotary  component  phase  angles  rotate 
CCW  with  depth  but  CW  rotary  component  phase  angles 


rotate  CW  with  depth  and  therefore  (3)  the  rotating  veetors 
will  align  at  nearly  the  same  orientation  but  the  timing  of 
the  alignment  (phase)  will  advanee.  In  short,  differing  decay 
scales  for  rotary  component  amplitudes  cause  effect  1,  and 
rotation  of  rotary  component  phases  either  together  or 
oppositely  eause  effect  2  or  3. 

[6i ]  Both  constant  Az  theory  and  linear  A:  theory  repro¬ 
duce  such  depth  changes.  However,  in  general,  constant  Az 
theory  does  not  mateh  well  the  observed  curvature  of  the 
semi  -major  axis  and  semi  minor  axis  profiles,  and  linear 
Az  theory  (with  a  constant  k  mueh  larger  than  molecular 
viscosity)  matches  better  the  observed  structure  of  all  the 
most  important  tidal  parameter  changes  in  the  Ekman  layer. 
The  average  scale  height  (where  fiz  =  k)  for  optimized  f3 
and  k  values  found  in  this  study  was  2.9  m.  This  implies 
an  extended  region  away  from  the  bottom  where  eddy 
viscosity  is  relatively  high  and  approximately  constant.  Of 
course,  since  the  ADCP  measurements  do  not  extend  down 
into  this  layer,  the  aetual  structure  near  the  sea  bed  eould  be 
different  and  Az  eould  have  a  different  form  than  the  one 
that  fits  tidal  velocities  2.5  m  and  higher  off  the  bottom. 
The  VR4  mooring  was  configured  with  a  higher-frequency 
ADCP  and  it  measured  velocities  as  elose  as  1  m  off  the 
bottom  as  part  of  its  bottom  depth  cell,  yet  error  levels  were 
not  anomalous  (Figure  5d)  compared  to  other  sites  with 
measurements  farther  off  the  bottom. 
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Figure  16.  Vertical  structure  of  various  estimates  of  Az  at  sites  SS5,  SS8,  CP3,  and  VR6.  The  theoretical 
linear  fits  to  observations  arc  shown  for  the  Kj  (orange),  M2  (black),  and  S2  (magenta)  constituents.  The 
time  mean  values  for  NCOM  (red)  and  Quoddy  (blue)  are  drawn  with  dots  to  indicate  the  depth  levels 
where  Az  is  defined.  The  estimated  gains  at  the  Kj  (green)  and  M2  (cyan)  tidal  frequencies  are  drawn  as 
solid  lines  for  NCOM  and  as  dashed  lines  for  Quoddy.  Gains  with  estimated  errors  that  are  greater  or 
equal  to  50%  of  the  gain  value  arc  not  drawn. 


[62]  Bottom  roughness  is  one  possible  explanation  for 
high  eddy  viscosity  near  the  bottom.  If  the  roughness 
elements  extend  beyond  a  hypothetical  viscous  sublayer, 
molecular  viscosity  is  unimportant  even  near  the  bed 
[ Soulsby ,  1983].  Of  course,  the  moorings  themselves  could 
act  as  roughness  elements.  As  described  by  Perkins  et  al. 
[2000],  BARNY  mounts  are  shaped  like  large  barnacles 
with  a  2  m  diameter  circular  footprint  and  0.5  m  maximum 
height  The  effect  of  a  BARNY  mount  on  the  tidal  currents 
is  unknown,  but  some  perturbation  of  the  flow  field  near  the 
sea  bed  should  be  expected. 

[63]  Werner  et  al.  [2003b]  used  Benthic  Acoustic  Stress 
Sensor  tripod  deployments  at  76  m  depth  on  the  southern 
flank  of  Georges  Bank  to  determine  M2  tidal  friction  shear 
velocities  (w*)  through  logarithmic  fits  of  observed  near¬ 
bottom  velocities.  They  used  standard  logarithmic  layer 
theory  to  derive  bottom  stress,  a  bottom  drag  coefficient, 
and  bottom  roughness.  M2  bottom  shear  stresses  estimated 
from  their  experiment  averaged  1.6-10  4  m2/s2  [Werner  et 
al ,  2003a],  1 1  times  larger  than  the  average  from  the 
northern  Adriatic  best  fit  linear  Az  theory  results.  This  factor 
is  similar  to  the  average  squared  ratio,  15,  of  M2  current 
semi  major  axis  amplitudes  in  their  study  and  this  one,  so, 
as  expected,  w*  scales  roughly  as  the  current  speed.  Inserting 


our  time  mean  w*  values  in  equation  (C2)  from  Werner  et  al 
[2003a],  gives  a  logarithmic  layer  thickness  of  ~2.0  m, 
which  is  below  the  depth  range  observed  by  most  ADCPs 
used  in  our  study.  For  unstratified  conditions,  Werner  et  al 
[2003a]  found  good  agreement  between  measured  tides 
above  the  logarithmic  layer  and  a  constant  Az  model  (linear 
Az  below  the  log  layer).  But  their  model  did  not  produce  a 
distinct  bulge  in  semi  major  axes  and  semi  minor  axes  as 
the  Adriatic  fits  did,  probably  because  they  used  Az  values 
(~0.03  m2/s)  approximately  40  times  the  northern  Adriatic 
best  fit  values,  with  an  associated  Ekman  depth  >100  m 
spreading  the  bulge  over  a  much  larger  depth  range 
However,  the  use  of  such  large  constant  Az  values  produces 
poor  fits  to  the  observed  Adriatic  tides. 

[64]  When  upper  water  column  biases  arc  removed  from 
NCOM  and  Quoddy,  (particularly  the  amplitude  bias  in 
Quoddy  and  the  phase  bias  in  NCOM)  they  too  match  the 
observed  currents  in  the  bottom  Ekman  layer  well.  The 
Mellor-Yamada  level  2.0  and  2.5  turbulence  closure  schemes 
implemented  in  these  models  seem  to  reproduce  vertical 
eddy  coefficients  that  are  well  suited  for  mimicking  realistic 
tidal  changes  near  the  bottom.  Although  the  errors  for  use 
of  linear  Az  theory  are  smaller,  it  is  an  optimized  best  fit 
from  site  to  site  and  constituent  to  constituent  which  is 
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Figure  17.  Cumulative  sum  versus  frequency  for  the  convolution  integral  between  the  Discrete  Fourier 
Transform  of  Az  and  the  Discrete  Fourier  Transform  of  at  sites  SS4  and  SS9.  The  colors  indicate 
NCOM  M2  (red),  NCOM  Ki  (cyan),  Quoddy  M2  (magenta),  and  Quoddy  K\  (blue).  The  solid  lines  are 
from  the  level  closest  to  1.5  m  off  the  bottom  (circles  mark  the  zero  frequency  values),  and  the  dashed 
lines  are  from  the  level  closest  to  15  m  off  the  bottom  (crosses  mark  the  zero  frequency  values).  In  each 
ease,  the  component  of  7\  with  the  strongest  stress  is  the  one  shown. 


reasonable  because  of  the  potential  for  spatial  changes  in 
bottom  roughness,  but  is  a  luxury  not  enjoyed  by  NCOM  or 
Quoddy.  However  the  fact  that  optimized  linear  Az  theory 
consistently  gives  much  better  results,  specifically  for  semi¬ 
diurnal  phase  and  diurnal  ellipse  orientation,  suggests  a 
dynamical  difference  Linear  Az  theory,  with  use  of  a  rela¬ 
tively  high  value  for  k ,  seems  to  be  particularly  adept  at 
matching  such  angular  changes.  The  sizes  of  the  theoretical 
bottom  shear  stress  ellipses  in  Figure  18  are  primarily 
dictated  by  the  values  of  k  and  if  these  extrapolated  results 
are  valid,  they  suggest  that  the  bottom  shear  stresses  of  the 
two  models  are  too  low. 

[65]  One  possible  explanation  is  that  use  of  equations  (1) 
and  (2)  assumes  a  logarithmic  velocity  profile.  Clearly 
equations  (6)  (11)  depart  from  this  in  agreement  with  the 
general  finding  of  Soulsby  and  Dyer  [1981],  who  suggest 
that  in  an  accelerating  (but  nonrotating)  tidal  flow  the  near- 
bed  velocity  profile  departs  from  the  usual  logarithmic 
form.  In  contrast,  Lueck  and  Lu  [1997]  in  their  tidal  channel 
measurements  found  that  the  velocity  profile  above  3  m  off 
the  bottom  in  the  along  channel  direction  matched  a 
logarithmic  form  well  and  that  departures  from  their  fits 
were  not  consistent  in  this  depth  zone  with  Soulsby  and  Dyer 
[1981  ]  acceleration  corrections.  But,  they  also  found  that  the 


velocity  profiles  in  the  across  channel  direction  matched  a 
linear  profile  instead  of  a  logarithmic  one.  Figure  18  from 
our  results  shows  that  bottom  stress  estimates  from  the 
quadratic  law  are  particularly  low  relative  to  estimates  from 
equation  (12)  in  the  semi  minor  ellipse  direction.  Another 
effect  to  consider  is  form  drag  [Chriss  and  Caldwell ,  1982] 
which  Lueck  and  Lu  [  1 997],  Ullman  and  Wilson  [  1 998],  and 
Werner  et  al.  [2003b]  all  used  to  explain  their  findings  of 
larger  roughness  lengths  or  drag  coefficients  than  typical 
values.  However,  it  is  less  clear  how  this  could  explain  our 
general  finding  of  implied  higher  bottom  stress,  as  the  sites 
near  Italy  had  very  muddy  bottoms  and  typical  bed  forms 
responsible  for  form  drag  are  not  expected  to  be  present. 

[66]  It  is  interesting  to  note  that  despite  a  wide  spread  in 
Az  gradients  between  the  various  model  and  theoretical 
solutions  shown  in  Figure  16,  all  the  solutions  seem  to 
converge  to  close  to  the  same  Az  values  in  the  bottom  1  m.  If 
the  vertical  derivative  is  carried  through  the  friction  term  in 
equation  (3),  then,  for  a  given  depth  level,  it  is  the  numerical 
value  of  A ,  that  multiplies  and  the  slope  of  Az  that 
multiplies  which  together  with  the  convergence  of 
solutions  to  a  particular  Az  value,  suggests  the  relative 
importance  of^r  for  determining  the  tidal  structure  in  this 
region.  If  the  velocity  structure  above  this  level  is  roughly 
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Figure  1 8.  Kinematic  bottom  shear  stress  from  theory  and  models  at  sites  SS8  and  VR6  for  (left)  M2  tides 
and  (right)  Ki  tides.  Red  is  from  Quoddy,  green  is  from  NCOM,  blue  is  from  a  model  with  tuned  linearized 
friction,  and  black  is  from  the  best  fit  linear  Az  theory.  The  black  solid  lines  are  from  equation  (12),  while 
the  black  dashed  lines  use  quadratic  drag  law  from  the  1-m  height  theoretical  results  with  Cd  chosen  to 
yield  the  same  maximum  rjp  magnitude  as  the  solid  ellipses.  Dots  indicate  relative  phasing.  The  direction 
of  rotation  of  Tolp  with  time  is  from  the  dots  to  the  gaps. 


logarithmic,  then  and  should  have  opposite  signs,  and 
therefore  higher  numerical  values  of  Az  in  Quoddy  and 
NCOM  could  offset  the  effects  of  larger /C  slopes  and  thus 
explain  how  the  diverging  Az  depth  profiles  result  in  similar 
tidal  solutions. 

[67]  The  analysis  of  the  time  structure  of  Az  from  Quoddy 
and  NCOM  and  the  fact  that  theoretical  fits  with  time- 
constant  Az  match  the  observations  best  show  that  time 
variation  of  Az  is  not  very  important  for  Adriatic  tidal 
momentum  balances.  As  the  CW  rotary  component  of  the 
Adriatic  tidal  flows  decays  slower  toward  the  bottom  than 
the  mean  flow,  tidal  components  become  more  dominant  in 
the  spectra  of  Az  and  and  the  interaction  of  the  time  mean 
of  Az  and  the  tidal  components  of  r.£  predominately  deter¬ 
mine  the  frictional  tidal  shear  stresses  in  this  region  where 
the  stresses  tend  to  have  maximums.  For  models  without 
mean  or  low-frequency  flows,  such  as  this  application  of 
Quoddy,  the  role  of  the  time  mean  of  Az  is  even  greater  in 
determining  the  frictional  shear  stresses  for  tides. 

[68]  Shear  stress  in  the  eross-axis  direction  of  the  Adriatic 
is  particularly  important  for  diurnal  tides.  Diurnal  tidal 
currents  near  the  bottom  rotate  with  depth  more  toward  this 
direction  than  semidiurnal  tidal  currents,  and  the  diurnal 


tidal  waves  propagate  in  this  direction  from  the  northeast  to 
the  southwest  eoasts  unlike  the  semidiurnal  tide  waves, 
whieh  propagate  along  the  axis  of  the  sea  [Malacic  et  a/., 
2000].  These  faets  may  aeeount  for  the  matching  bias 
structure  for  K!  amplitudes  in  NCOM  and  Quoddy.  Both 
models  use  spatially  uniform  bottom  roughness  and  there¬ 
fore  may  miss  the  potentially  significant  effect  on  Kj  of  a 
northeast  to  southwest  bottom  drag  difference  caused  by 
varying  sediment  types  (sand  to  mud).  The  trend  for  the  SS 
line  in  Figure  8  (top  left)  is  in  the  right  sense  for  this  where 
a  moderate  and  uniform  model  Cd  would  not  damp  enough 
energy  near  the  northeast  coast  but  damp  too  much  near  the 
southwest  eoast.  Lower  bottom  shear  stress  in  the  semi 
minor  axis  direction  may  also  have  an  effect. 

7.  Conclusions 

[69]  A  large  observational  database  unmatched  in  previous 
studies  was  compiled  and  used  with  a  suite  of  mathematical 
and  numerical  models  in  an  effort  to  explore  Adriatic  tidal 
dynamics.  Fifteen  bottom-mounted  ADCPs  deployed  for 
more  than  six  months  in  the  northern  Adriatic  during  a  time 
period  with  generally  little  stratification  were  able  to  resolve 
strong  tidal  current  structure  changes  at  1-m  intervals  from 
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heights  of  3  m  off  the  bottom  upward  to  heights  where  the 
tidal  current  structure  was  nearly  vertically  uniform.  Tidal 
currents  in  much  of  the  water  column  were  near-reversing 
tides,  with  M2  currents  having  the  largest  amplitudes  and  S2 
and  K,  currents  both  having  significant  energy.  Near  the 
bottom,  tidal  current  ellipses  were  all  shortened  and  broad¬ 
ened,  semidiurnal  currents  led  upper  water  column  currents, 
and  diumal  tidal  current  ellipse  orientations  rotated  CCW 
descending  through  the  boundary  layer. 

[70]  Such  changes  are  in  accord  with  tidal  bottom  Ekman 
layer  theory  and  match  solutions  using  a  linear  form  for  Az 
that  are  obtained  here  for  conditions  of  no  slip  at  the  bottom 
and  no  stress  at  00.  Linear  Az  solutions  have  30%  smaller 
errors  than  constant  Az  solutions,  but  only  if  the  constant 
term  in  the  linear  form  is  far  above  molecular  values  NCOM 
and  Quoddy  simulations  were  also  compared  to  observations. 
Accounting  for  upper  water  column  biases,  the  model 
simulations  in  the  bottom  Ekman  layer  accurately  reproduce 
the  major  tidal  structure  changes  with  depth,  although  not 
as  well  as  optimized  theory. 

[71]  Analysis  of  the  time  structure  of  Az  solved  by  using 
Mellor-Yamada  turbulence  closure  in  NCOM  and  Quoddy 
show  that  the  time  mean  of  Az  predominately  controls  the 
frictional  tidal  shear  stresses.  Comparison  of  the  depth 
structure  of  Az  between  the  model  solutions  and  the  theo¬ 
retical  fits  shows  a  large  spread  of  values  in  the  bottom  24-m 
layer,  but  also  shows  convergence  toward  similar  values  at 
depths  less  than  1  m  off  the  bottom.  The  optimized  theory 
fits  show  that  the  constant  part  of  Az  was  most  important  in 
matching  the  tidal  orientation  and  phase  changes  and  a 
nonzero  Az  slope  was  important  in  matching  the  tidal 
amplitude  changes. 

[72]  A  comparison  of  kinematic  bottom  shear  stresses 
from  the  optimized  theory,  from  Quoddy,  from  NCOM,  and 
from  a  data  assimilation  model  using  an  optimized  linear 
friction  parameter  shows  a  wide  range  of  values  produced 
by  these  different  models  in  simulating  the  northern  Adriatic 
tides.  Quoddy  tends  to  have  the  smallest  stresses,  but  NCOM 
stresses  arc  similar.  Theoretical  stresses  are  larger,  especially 
in  the  cross-axis  direction.  The  model  using  linearized 
friction  has  the  largest  stresses  in  the  along-axis  direction 
but  weakest  in  the  cross-axis  direction.  For  near- reversing 
tidal  currents,  using  a  quadratic  drag  law  or  a  linear  friction 
parameter  can  produce  low  bottom  shear  stresses  in  the 
semi  minor  axis  direction. 

[73]  Further  work  is  needed  to  explore  the  impact  of 
potentially  underestimating  this  direetional  component  of 
stress  in  numerical  models  of  the  Adriatic  and  elsewhere. 
Also,  the  northern  Adriatic  with  its  general  lack  of  stratifi¬ 
cation  during  half  of  the  year  and  significant  semidiurnal 
and  diumal  tides  would  be  a  good  location  to  test  the  validity 
of  equations  (6)  (11)  near  the  bed  from  examination  of 
appropriate  current  measurements  from  long-term,  bottom 
tripod  deployments. 
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